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^ ! Abstract 

• Real time thermalization and relaxation phenomena are studied in the low energy 

density phase of the 2+1 dimensional classical 0(2) symmetric scalar theory by 
solving numerically its dynamics. The near-equilibrium decay rate of on-shell waves 
and the power law governing the large time asymptotics of the off-shell relaxation 
agree with the analytic results based on linear response theory. The realisation of 
the Mermin- Wagner theorem is also studied in the final equilibrium ensemble. 

^ : 

The approach to equilibrium of initially out-of-equilibrium states is a highly important 
Q-i! issue in many branches of physics ranging from inflationary cosmology (the spectrum 

of density fluctuation in the early universe) through particle physics (the problem of 
bariogenesis, formation of DCC in heavy ion collisions) to statistical physics (dynamics 
of phase transitions, realization of Boltzmann's conjecture []. 

In the recent field theoretical literature thermalization and relaxation of classical fields 
is intensively investigated Its interest follows from the presence of bosonic degrees of 

freedom with high occupation numbers, for instance, in cosmological applications 0. In 
fact, efficient methods for the determination of the exact, real time evolution, e.g. numer- 
ical integration of the equations of motion are known only for classical fields. Truncation 
and expansion schemes should be benchmarked against this exact solution. 

Recently, evolution of equal-time 1PI correlation functions derived from the effective 
time dependent action [[| were confronted with the results of the exact time evolution. 
By solving equations non-local in time, obtained from 2PI effective action thermalization 
of quantum fields was demonstrated ||. These approximate methods can be formulated 
both for classical and quantum cases. In the quantum case the analogue of the classical 
ensemble averaging over the initial conditions is the quantum expectation value. With 
this correspondence the formal derivations, and hence the results, are quite similar ||. 

*mazsx@cleopatra. elte.hu 
tszepzs@antonius.elte.hu 

1 An ensemble of isolated interacting systems approaches thermal equilibrium at large times. 
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In statistical physics real time studies of complex nonlinear systems have shown that 
expectation values of observables are correctly reproduced by averaging over the micro- 
canonical evolution of a single equilibrium system, even if a very small part of the phase 
space is covered by the motion during the interval of observation [jj . 

In this Letter we wish to gain insight into the validity of the linear response theory, by 
comparing its results to exact solutions. A similar attempt was already done in ||. There 
the early, far from equilibrium time evolution was confronted with the linear response 
results, and a relaxation slower than expected has been observed. We will show later 
(see Fig. that the decay in the linear regime realizes the fastest relaxation rate of this 
quantity during the whole time evolution. Exploration of the range of validity of the 
linear regime is possible only by following the evolution very long, and averaging over 
initial conditions in order to raise the signal over the noise level. 

The model and its numerical solution. We concentrate on the 0(2) symmetric classical 
scalar theory with a small explicit breaking term in its Lagrangian: 

£ = + \(d^) 2 - \m^\ - \m^\ + + M>x, (1) 

where h controls the explicit symmetry breaking, and m 2 < 0. We have studied the time 
evolution of a system, discretized on lattice in two space dimensions at so low energy 
density (see below) that according to the mean field analysis this would correspond to the 
broken symmetry phase. 

The dynamics eventually drives the system towards thermal equilibrium. For h ^ 0, 
the equilibrium state has large magnetisation, nearly corresponding to the minimum of the 
classical potential. Fluctuations around this state are naturally divided into Goldstone 
(light) and Higgs (heavy) excitations. These excitations experience an effective potential, 
which agrees very well with the result of finite volume perturbation theory. The relaxation 
into this state is the main subject of this paper. We shall present a detailed comparison 
of the exact time evolution with the linear response theory. 

A second relaxation process can be initiated from the equilibrium with h ^ 0, if the 
magnetic field is switched off. The final h — equilibrium, however, obeys in the thermo- 
dynamical limit the Mermin-Wagner-Hohenberg theorem || , which states the absence of 
spontaneous magnetisation in two dimensional systems with continuous symmetry. For 
equilibrium two dimensional finite volume systems even for vanishing external source, at 
very low temperature there is a non-zero magnetisation with a well defined direction as 
shown in ||. The finite volume magnetisation selects the direction with respect to which 
we can define a parallel and a perpendicular mode. This natural choice of base suggests 
the use of the Goldstone boson terminology in the h = case too. The way the finite 
magnetisation disappears as the volume goes to infinity will be touched upon shortly at 
the end of our discussion. 

The exact time evolution of the lattice system was studied by introducing dimension- 
less field variables <3>j — > y/a&i and rescaling further these fields like $rescaled — $ 12 -\/A/6 

and h Tesca ^ = /ia/A/6. This sets the bare coupling to A resca ^ e< ^ = 6. We have chosen 
m 2 a 2 = —1, therefore all dimensionful quantities are actually expressed in units of \m\. 
Squared lattices of size ranging from 64 x 64 to 512 x 512 were used. 

The explicit symmetry breaking parameter h is chosen in the range ... 0.0025/^6. 
Initially the system is at rest near the origin of the <3>-space. A small amplitude white 
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Figure 1: Time evolution of |k| power spectra of Higgs (left) and Goldstone (right) field 
components (512 x 512 lattice, h = 0.0025/^ ) 



noise determines the field in every lattice point ($i(x, t = 0) = ?7i(x), $2(x, t = 0) = 
r] 2 (x), $!(x, t = 0) = $ 2 (x,t = 0) = 0, (r]i(-x)r]j(y)) ~ 5 x ,y^,j)- This means, that the 
space average of the two-component order parameter (OP) ($i , $2^) is very close 
initially to zero, which is (approximately) a local maximum of the effective potential. 
The roll-down towards the minimum is oriented by the external field h. It is the potential 
energy difference of the initial and the final states, which is redistributed between the 
modes. The magnitude of the noise was so small, that the final energy density of the 
system is exclusively determined by this difference in the potential energy. 

Numerical diagonalization of the measured 2x2 matrix of equal time fluctuations 

V V V 

(Cij(t) = $j(x, £)<&j(x, t) — $j(x, t) $j(x, t) ) gave as uncorrected degrees of freedom 
the radial (Higgs) and angular (Goldstone) components of ($ 1; $ 2 ) : 

$2(X,*) 



$ H (x, t) = \/<5> ?(x, t) + $|(x, *), $ G (x, t) = $ H ,o arctan ^ f y (2) 



§ H0 = \m\y/6/X + h/2\m\ 2 + 0{h 2 ) being the minimum of the bare potential valley. 



-V rrV- ^-rV 



This statement means that & H (t)$ G (t) $ H (t) ®c(t) for all t. The boundedness of 
$2 fluctuations follows from the IR-regulating effect of the explicit symmetry breaking 
term. The position of the tree level minimum scaled by |m| is &h,o = 1 + h/2 + 0(h 2 ), 
and the scaled squared masses read as m 2 H = 2 + h/2 + 0{h 2 ) for the radial mode and 
m G = h + 0(h 2 ) for the angular one. 

Relaxation to the h ^ equilibrium. During and immediately after the roll-down the 
main mechanism of the fast excitation of spatial fluctuations is the parametric resonance, 
as observed and studied by many authors, e.g. |I0, DJ. Also the dephasing of different 
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oscillation modes can be observed as reported in [|, [L2|. The time scales of this stage 
are characterised by the cutoff and by the masses of the radial and angular modes. The 
measurement of these latter is discussed later. 

The early time evolution can be characterised by the variation of the kinetic power 
spectra of the Higgs and Goldstone fields (see Fig. p. After dephasing the potential and 
kinetic energy of each mode is balanced in itself, therefore the energy content of different 
modes may be characterised by the kinetic power spectrum 

n k 

P H , G (\k\) = &„ Q , (3) 
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where the averaging should be taken over the polar angle in momentum space. In the 

final equilibrium state this quantity does not depend neither on |k|, nor on which field it 

refers to. Its value corresponds to the temperature, which was measured in our case to 

be P h ,g(M) = Tkin = °- 125 =•= 0.0005. This kinetic temperature can obviously be defined 



also out-of-equilibrium as T kin = <& 2 H G . 

The initial parametric peak in the power spectrum of the Higgs field disappears rapidly. 
It proceeds mainly through decays into pairs of the lighter Goldstone modes. This damp- 
ing is reflected by an energy transfer from the Higgs to the Goldstone modes. In addition, 
the excitation of spinodally unstable Higgs modes also contributes to the relaxation, as it 
was shown in ||13[ . 

The parametric low |k| Goldstone peak, however, survives as long as t = 10 4 ...10 5 . 
During its long relaxation the energy is slowly transferred back to the Higgs modes. 
This process qualitatively corresponds to off-shell emission of Higgs waves. The final 
equilibrium is reached in a slow, non-exponential relaxation of Tq u to T^. 

The mechanism sketched above is suggested by the structure of the analytic formu- 
lae derived in the perturbative relaxation analysis, which should be relevant at least to 
the large time asymptotics of the evolution. The forward time evolution of an initial 
configuration is determined by the classical self-energy function n(x, t): 
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(4) 



where z(k) = ik F(k) — -P(k), determined by the corresponding initial configuration 
F(k) = $(t = 0, k), P(k) = dt$(t = 0, k). 

For the calculation of the Fourier transform of n(x, t) the procedure described in 
ni [T5I] was used. (For the quantum treatments, which lead for large occupation numbers 



to the same results, see |16|, |P7|1). 

First one makes the shift $i — > $i + $ in order to describe the constant equilibrium 
background, where $ is the classical ensemble average of $#(x). Its direction is actually 
selected by the explicit symmetry breaking. Then in the linear approximation one finds 
$ff(x) ~ $i(x) and $g( x ) ~ < ^2( x ) < ^ ) h,o/5 ) - We refer to our previous publication 
Appendix B] for the detailed derivation of the self energies, which read now for the Higgs 
and Goldstone fields, as follows: 
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where uj 2 H/G k 



m 



H/G 



|k| 2 . The relevant mass values were calculated numerically from 
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the relaxed field configurations with the method developed in [13]. At the parameter 
values specified in our simulations the typical masses are am H rs 1.3 and am G < 0.033. 

The imaginary part of the self energy, which accounts for the damping phenomena, may 
be extracted using the principal value theorem. There is a technically important difference 
between our present formulae and the ones derived in |Tj|. It follows from the fact that 
our system is defined in two space dimensions. This circumstance modifies the volume 
element in the momentum space and an integrand of more complicated analytical structure 
appears when calculating the imaginary part of the classical self energy function. 

The explicit expressions of the imaginary parts are the following: 
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It is the off-shell damping coming from the contributions of the cuts to (^) which 
has a direct impact on the late time evolution of the OP we are analysing. We followed 
the method, developed in 
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which determines the leading power law tail of the 
relaxation. Contrary to the 3 + 1 dimensional case in 2 + 1 dimensions the integrand 
of Eq. (H) has extra cuts. These can, however, be directed in a way that either they do 
not contribute (Goldstone case) or their contribution is suppressed by a factor of e~ mHt 
(Higgs case). The large time asymptotics of the two field components is the following: 
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Figure 2: Relaxation of Goldstone spatial fluctuations with time on (e-based) log-log scale 
plot. (256 x 256 lattice, h = O.OOOl/v^, m G = (62 ± 5) • 10" 4 , m H = 1.329 ± 0.001 • lO" 3 , 
the average over 69 runs is shown.) 
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where the threshold frequency values, appearing above are given by Q 2 H / G k = k 2 + rn? H i G 
and ^l 2 H±GM = k 2 + (m H ± m G ) 2 . 

In order to obtain numerical evidence for the presence of the predicted power law tail, 
we measured the quadratic spatial fluctuation moments of both fields as a function of 

time. They are defined as Fluct(t) = $(x, t) 2 — ( $(x, t) \ . Substituting (|Tl| ) into this 

definition we get ~ t~ 2 like relaxation for late times. This behaviour is indeed observed 
as the fitted power shows in Fig. for the example of the Goldstone mode. 

The formulation of the theory upon which the above formulae were obtained actually 
uses ensemble averages over the initial data of the classical evolution. We found, that 
although every single run seems to relax even quantitatively in the same way, this late 
time evolution can be extracted from the noisy data only, if this averaging is performed 
indeed (116 runs were involved). We could recognise a ~ 1/t 2 decay for Higgs modes 
too, but because of the low signal/noise rate we could not do quantitative analysis. In 
particular the analysis is made difficult by the fact, that Higgs modes decay rapidly and 
vanish in the noise before the power law tail would be reached. 

Relaxation behaviour of OP (<E># (k=0)) was also investigated, and — in accordance 
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with the expectation above — an oscillation damped by ~ £-i±o.o2 was f oun d. Moreover, 
this oscillation, is observed around a value, monotonically approaching its equilibrium 
value as ~ t~ 2±0 01 . This latter behaviour is explained by the fact, that the exact equation 
of motion for OP contains Fluct(t) as a time dependent parameter |10, Its ~ t~ 2 



like behaviour is inherited by the slowly varying part of OP. (The errors of the exponents 
come from averaging over 92 runs on a 256 x 256 lattice.) 

The h equilibrium. When one follows the evolution of the system for long times 
(t = 10 6 ), there seems no further relaxation to take place, i.e. the initial signal has been 
lost in the thermal noise. In order to convince ourselves of having arrived to thermal 
equilibrium we compare the measured values of , and the masses mo and m# to the 
analytical estimates of these quantities coming from the one-loop lattice effective potential 
evaluated for the same size lattice as used in the simulation. In the analytic expression 
we used the measured kinetic temperature. The equilibrium masses mg and m# were 
determined in the simulation both from the corresponding correlation functions and by 
fitting the oscillatory motion around the equilibrium of the corresponding OP-components 
as described in fll3| . 

We have checked that the system with the present initial conditions is deeply in the 
Coulomb phase i.e. the temperature was about four times smaller than the Kosterlitz- 
Thouless critical temperature T^t- in the absence of explicit symmetry breaking on 
the critical line between T = and Tkt the Goldstone correlation length is expected 
to diverge with the lattice size. The explicit symmetry breaking parameter h acts as 
an infrared regulator. Indeed, the measured Goldstone correlation length is found to 
be proportional to L, but for the largest sizes (L = 256), where the IR cutoff h begins 
to dominate. On the other hand, the inverse correlation length in the Higgs direction 
is finite, its value being 3% smaller than the two-loop mass and 6% smaller than the 
one-loop value. 

A very good agreement of with $ is found where $ is given by the minimum of 
the effective potential in the radial direction, the relative deviation being O(10~ A ). The 



discrepancy between the measured (with the method described in [13 ) and calculated 
masses was less than 1% and 5% for the Higgs and Goldstone modes, respectively, the 
perturbative values being systematically smaller. Our measurements of the Goldstone 
mass became very noisy for small values of h < 0.00001/ V6. 

In the equilibrium system one can proceed to "experiments" , which check the correct- 
ness of the on-shell decay rates computed in the linear response theory. These are the 
simple zeros of the denominator of Eq.([|) which determine the exponential damping of 
on-shell excitations. The damping rates are obtained by substituting |k| = /c into Eqs. 
(0) and (H). To leading order in (A$/6) 2 T the rates read (assuming m H > 2m G ) as 

IuiII^ 2 = m|) _ /A<5\ 2 T 



= ^- ^ = — - == (12) 
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(13) 



On equilibrium configurations single k-modes have been superimposed with an ampli- 
tude between 0.01 . . . 0.2. Time evolution of these excited modes showed perfect exponen- 
tial decay. The exponents did, however, depend on the amplitude. For small amplitudes 
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Figure 3: Numerical and analytical values for Goldstone on-shell decay rate as a function 
of |k| (128 x 128 lattice, h = 0.0025/^, m G = 0.0319, m H = 1.33, $ = 0.96, T = 0.125) 



the fit was unstable, the mode was quickly lost in the noise. For bigger values nonlinear 
effects did show up. Therefore we were looking for a plateau in the amplitude dependence 
of the decay rates between these two extremes. Then its value has been compared to the 
analytical estimates (see Fig .[3]). The error bars show how well-defined plateau could be 
found. This numerical experiment was performed both for Higgs and Goldstone waves, 
but the errors for Higgs decay were too large again. The value of Th computed from 
the linear response theory was reproduced within an accuracy of ±50%, but we can not 
say anything about its predicted k-dependence. The reason for this seems to be, that 
the evolution becomes very soon nonlinear as we try to increase the excitation amplitude 
higher that the noise level. 



As was indicated in [21| finite size effects may be highly important in systems with 
Goldstone-like excitations. This circumstance made necessary the use of finite volume 
perturbation theory in the analysis above. Neither the lattice summed perturbative effec- 
tive potential, nor the numerically found mass values showed any L dependence if L > 64. 
We have checked for the L independence in all cases where the exact evolution and linear 
analysis were confronted (e.g. the values of damping rates). 

Diffusion to the h = equilibrium. The most care must be taken when the macroscopic 
magnetisation M(L) is considered. For finite volume, the ensemble average computed 
from the low temperature spin-wave approximation is non-zero || even in the absence of 
explicit symmetry breaking: 

M\L) = ($/) 2 + ($f ) 2 « (2L)- T ^ 2 \ (14) 

For each realization of the canonical ensemble observed in Monte Carlo simulations the 
magnetisation vector has a well-defined direction 9 = arctan ( $ 2 /$! 



It is an interesting question, by what mechanism the MWH-theorem is realized as 
L — > oo. It was shown in Monte Carlo simulations that the erasure of the order is realized 
by the diffusion like displacement of the direction of 9 ||. 

In our real time study, for the investigation of the onset of the finite volume version 
of MWH-theorem a number of runs were continued from the equilibrium state reached 
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with h 7^ 0, after the magnetic field was switched off. In each individual system OP 
begins to circle around the origin with an average radius M(L). This motion is probably 
due to a random nonzero angular momentum of the initial configuration. Subtracting 
the angular momentum of OP, the relevant diffusion-like motion remains. MC studies 
employing first order Monte Carlo "time" -evolution of non-equilibrium one dimensional 
systems have shown that the sign for SSB is the exponential decrease of the diffusion 



constant with L \T2\. In accordance with the expectations based on the MWH-theorem 
we find in the present case that the diffusion constant decreases according to a power law 
with the exponent —1.16 ± 0.1. 

In conclusion we can state that the numerical study of the Hamiltonian dynamics 
of the classical 0(2) symmetric scalar model in 2+1 dimensions provides a non-trivial 
check of our understanding of real time relaxation phenomena. Numerical results both 
for the late-time OP-asymptotics and for the decay rate of on-shell waves (with well- 
defined wave vector k) were found to be in agreement with the linear response theory. 
In the h equilibrium the agreement of the masses coming from the perturbatively 
calculated effective potential with the numerically established excitation masses was also 
verified. Finally, we have demonstrated the real time manifestation of the finite volume 
MWH-theorem by measuring the large L asymptotics of the angular diffusion rate of the 
macroscopic order parameter. 
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